%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% File:               least_sq.m
%
% Authors:            Sergio Ascencio and Miguel Rueda
%
% Description:        Computes objective function (LS Pesendorfer's approach)
%
% Language:           MATLAB R2013b (8.2.0.701) 64 Bit
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


function [diff,J]=least_sq(betas,X,phats)

%Extended
% varsx={'l_d_RC_others' 'l_margin' 'lcasillas' 'school_ic' 'l_diff_PRI_PAN' 'lpop' 'l_turnout' 'local' 'ldist_city' ...
% 'l_d_RC_others_sq' 'l_margin_sq' 'lcasillas_sq' 'school_ic_sq' 'l_diff_PRI_PAN_sq' 'lpop_sq' 'l_turnout_sq' 'ldist_city_sq' ...
% 'l_d_RC_others_l_margin' 'l_d_RC_others_lcasillas' 'l_d_RC_others_school_ic' 'l_d_RC_others_l_diff_PRI_PAN' 'l_d_RC_others_lpop' 'l_d_RC_others_l_turnout' 'l_d_RC_others_local' 'l_d_RC_others_ldist_city' ...
% 'l_margin_lcasillas' 'l_margin_school_ic' 'l_margin_l_diff_PRI_PAN' 'l_margin_lpop' 'l_margin_l_turnout' 'l_margin_local' 'l_margin_ldist_city' ...
% 'lcasillas_school_ic' 'lcasillas_l_diff_PRI_PAN' 'lcasillas_lpop' 'lcasillas_l_turnout' 'lcasillas_local' 'lcasillas_ldist_city' ...
% 'school_ic_l_diff_PRI_PAN' 'school_ic_lpop' 'school_ic_l_turnout' 'school_ic_local' 'school_ic_ldist_city' ...
% 'l_diff_PRI_PAN_lpop' 'l_diff_PRI_PAN_l_turnout' 'l_diff_PRI_PAN_local' 'l_diff_PRI_PAN_ldist_city' ...
% 'lpop_l_turnout' 'lpop_local' 'lpop_ldist_city' ...
% 'l_turnout_local' 'l_turnout_ldist_city' ...
% 'local_ldist_city' 'ldist_pri' 'gov_PRI' 'ldist_pan' 'gov_PAN'};

betas=reshape(betas,[length(betas)/4,4]);

p_PRI=phats{1,1}(:,1:2); %remember prob H in first, M in second column
p_PAN=phats{1,2}(:,1:2);

psi=best_response(betas,X,phats);

    Z_PAN=X(:,end-1:end);
    Z_PRI=X(:,end-3:end-2);
    X=X(:,1:9);

X_PRI=[X Z_PRI p_PAN ones(size(X,1),1)];
X_PAN=[X Z_PAN p_PRI ones(size(X,1),1)];

psi_PRI=psi{1,1}(:,1:2);
psi_PAN=psi{1,2}(:,1:2);

diff_PRI=p_PRI-psi_PRI;
diff_PAN=p_PAN-psi_PAN;

diff=[diff_PRI(:);diff_PAN(:)];

%Jacobian

auxPAN_M=[zeros(size(X_PRI)) zeros(size(X_PRI)) kron(ones(1,size(X_PAN,2)),-psi_PAN(:,2).*(1-psi_PAN(:,2))).*X_PAN  kron(ones(1,size(X_PAN,2)),psi_PAN(:,2).*psi_PAN(:,1)).*X_PAN];
auxPAN_H=[zeros(size(X_PRI)) zeros(size(X_PRI)) kron(ones(1,size(X_PAN,2)),psi_PAN(:,2).*psi_PAN(:,1)).*X_PAN  kron(ones(1,size(X_PAN,2)),-psi_PAN(:,1).*(1-psi_PAN(:,1))).*X_PAN];

auxPRI_M=[kron(ones(1,size(X_PRI,2)),-psi_PRI(:,2).*(1-psi_PRI(:,2))).*X_PRI kron(ones(1,size(X_PRI,2)),psi_PRI(:,2).*psi_PRI(:,1)).*X_PRI zeros(size(X_PAN)) zeros(size(X_PAN))];
auxPRI_H=[kron(ones(1,size(X_PRI,2)),psi_PRI(:,2).*psi_PRI(:,1)).*X_PRI kron(ones(1,size(X_PRI,2)),-psi_PRI(:,1).*(1-psi_PRI(:,1))).*X_PRI zeros(size(X_PAN)) zeros(size(X_PAN))];

J=[auxPRI_H;auxPRI_M;auxPAN_H;auxPAN_M];

end

